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Abstract 

A canard explosion is the dramatic change of period and amplitude of a limit cycle of a 
system of non-linear ODEs in a very narrow interval of the bifurcation parameter. It occurs 
in slow-fast systems and is well understood in singular perturbation problems where a small 
parameter epsilon defines the time scale separation. We present an iterative algorithm for the 
determination of the canard explosion point which can be applied for a general slow-fast system 
without an explicit small parameter. We also present assumptions under which the algorithm 
gives accurate estimates of the canard explosion point. Finally, we apply the algorithm to the 
van der Pol equations and a Templator model for a self-replicating system with no explicit small 
parameter and obtain very good agreement with results from numerical simulations. 


1 Introduction 

We consider singular perturbation problems of the form 

X = F{x,y,z,e), (1) 

y = G{x,y,z,e), 

where e is small parameter and 2 is a bifurcation parameter. All vector fields are assumed to be 
analytic in their arguments. In the singular limit e = 0 the system is assumed to have a critical 
manifold, that is, a manifold of critical points. We restrict our attention to planar systems with 
X G M and y G M such that the critical manifold is a curve. The parameter e defines a time-scale 
separation, and for singnlar perturbation problems in the special form 

X = F{x,y,z,e), (2) 

y = €H{x,y,z,e), 

the state variables x and y can be identified as the fast and the slow variable, respectively. For the 
system (2) the critical manifold is defined by F(x, y, z, 0) = 0, and it follows from Fenichel theory 
[7, 9] that, under certain regularity conditions, there exists a slow invariant manifold e-close to the 
critical manifold. For simplicity, we will suppress the explicit dependence of F and 77 on e in the 
system (2) from now on. 
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Canard explosion. The canard explosion is the dramatic change of amplitude and period of a 
limit cycle born in a Hopf bifurcation in a very narrow parameter interval. The phenomenon is well 
understood in singular perturbation problems of the form (2) [1, 17, 11, 3]. Canards are solutions 
which are e-close to a critical manifold on the form rjo = ?7o(y, z) defined from F(r/o(y, z),y,z) = 0 
both on an attracting part where dxF < 0 and a repelling part where dxF > 0. Thus, canards are 
slow manifolds with both an attracting and a repelling part. They only exist for a narrow range (of 
order of the parameter z, and the canard explosion occurs when a segment of a limit cycle 

is a canard. A unique asymptotic expansion of a parameter value z where the canard occurs, the 
canard point, and an expansion of the corresponding (maximal) canard can readily be obtained. 

However, canard explosions are also observed in slow-fast systems where there is no explicit 
small parameter that separates the timescales, 

x = F{x,y,z), (3) 

y = H{x,y,z), 

and where existing singular perturbation theory cannot be applied directly. An example is the 
Templator model which we consider in Section 3. In [5] a modihcation of the iterative method of 
Fraser and Roussel [8, 14], devised to construct regular slow manifolds, was proposed to determine 
a canard point for a general system of the form (3). The method was successfully applied to the 
van der Pol equations where it was shown that the first terms in the asymptotic expansion of the 
canard point are determined correctly, and on a Templator model for a self-replicating system with 
no explicit e where canard explosions have been found numerically. 

In this paper, we will consider a modified version of this method. Being based on a linearization, 
the method is explicit, in contrast to the original one. Furthermore, we present explicit non¬ 
degeneracy conditions as well as “smallness conditions” and show how they guarantee accurate 
approximations for the canard solution and the canard explosion point for the system (3). We 
illustrate the method with examples. 

Notation. All norms will be denoted by | • |. Superscripts with n G No will be used to denote 
partial sums such as; 

n 

= n>0, (4) 

i=0 

with each of the terms in the sum being enumerated through subscripts. Following this convention 
means that r]^ = r]o- We will in this paper suppose that the vector-field is analytic. In particular we 
will consider complex sets of the form x G (a, h) -|- iu and y G (c, d) + ia. Here we by (a, b) + 1^ C C 
and (c, d) + ia C C denote the complex v and cr-neighborhoods of the real intervals (a, b) resp. 
{c,d). We will then use \h\iy to denote the sup-norm of an analytic function h over the domain 
(a, b) + iv. This representation gives the following compact form of Cauchy’s estimate: 

Lemma 1 Let h = h{x) be analytic in x G (a, b) + in. Then 

\h'\u-^ < 0 < ^ < n. □ 


Proof See [15, Theorem 10.26]. ■ 

The iterative method of Fraser and Roussel. A slow manifold for Eqns. (2) of the form 
X = r]{y, z) fulfills the invariance equation 

-edyriH{y, y, z) + F{r], y, z) = 0, 


(5) 
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obtained by eliminating time. Fraser and Roussel [8, 14] developed an iterative method for the 
approximation of slow manifolds from the scheme 

-edyrf~^H{rf,y, z) + F{ri^, y, z) = 0, (6) 

starting from the critical manifold rp = r]Q{y,z). The function tj = satishes (5) up to the error 

-edyynH{y^,y,z), (7) 

using the notation (4). In [16] it was shown using Cauchy estimates to control (7) that this proce¬ 
dure, even for an arbitrary number of slow and fast variables, leads to slow manifolds exponentially 
close (0(e“'^/^)) to invariance. The equation (6) is non-linear in rp, but Neishtadt [13] showed that 
the same convergence is obtained if F and H in (6) are linearized with respect to x at such 
that the equation becomes linear in r]n = — y"'~^- 

A crucial assumption for the success of these methods is that 19x^1 e. This is violated 
near fold points of the critical manifold where dxF = 0, and this is where canards may occur. To 
study canards in systems with fold points one instead makes use of the fact that dyF p 0 and solve 
F{x, y,z) = 0 for y = Co(a^) z). A fold point (xq, Co(3;0) Mo)) is then characterized by 9xCo(a^O) Mo) = 0 
for some parameter value z = /Uq- For a slow manifold y = C{x, z) the invariance equation becomes 

-dxCF{xX,z) + eH{x,C,z) = (8) 

and following Fraser and Roussel it is proposed in [5] to solve this iteratively from 

-dxC~^F{x, C, z) + eH{x, C, z) = 0 (9) 

starting with the critical manifold = Co- A solution X{x, z) of (9) will generally have a singularity 
close to the fold, but this may be canceled by an appropriate choice oi z = yF. Then X is well- 
defined at the fold, and hence represents a canard. Thus, this procedure yields a sequence yF of 
approximations to the canard point as well as approximations y = CP{x,yF) of the corresponding 
canard. In the present paper we follow Neishtadt again, and linearize Eqn. (9) at such that 
it can be solved explicitly for CP. A precise formulation of this procedure requires some notation 
which we turn to now; the complete algorithm is described as Algorithm 1 below. 

2 The modified Fraser-Roussel algorithm for canards 

Consider the system (1) and let 


V{x,y,z,e) 


( F'ix,y,z,e)\ 
\G{x,y,z,e)J ’ 


denote the vector-held. A critical manifold for (1) is a smooth curve of hxed points V(x, y, z, 0) = 0 
within the (x,y)-plane for (1)^=0. We shall assume that this manifold can be parametrized by x so 
that y = Co(a^) z). We then have 

Lemma 2 The critical manifold is normally hyperbolic at (xq, Co(a^O) ^)) if and only if 


dyV • (-^xCo, 1) / 0. 


□ 
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Proof The layer equations of (1) are 

where x and y in general are both non-zero. An equilibrium (xq, Co(a^O) z)) is then normally hyper¬ 
bolic if the linearization, described by the Jacobian: 

DV = {d,v dyV), 

only has one zero eigenvalue. Since V{xXo{x,z),z,0) = 0 we have by implicit differentiation that 

d.v = -dyVd^Co- 

Therefore we can write DV as 

DV = {-dyVd^Co dyV ). 

The eigenvalnes of DV are then directly obtained; 

0 and dyV ■ {-dxCoA), 


and the resnlt therefore follows. ■ 

Following Lemma 2, the critical manifold y = Co{x,z) loses normal hyperbolicity at point 
(xo,Co(a^o,2:)) where 

dyV • (-5,Co, 1) = 0. (10) 

Generically (see (Al) below) snch a point is a fold point. 

Remark 1 Geometrically, (10) means that the critical fiber of (xq, Co(3;0)- 2 )) is tangent to the 
graph y = (o{x,z) at this point. □ 

One of the main aims of this paper is to present a method that applies to systems of the form 
(3) without an explicit e. To this end, we define a manifold y = Co{x, z) from 

F(x,Co(x,z),z) = 0. 

Were there an e multiplying ff, as in (2), this wonld be the critical manifold for e = 0 which can 
be nsed be used as a starting point for a Fraser-Roussel iteration as described in § 1. The precise 
assumptions needed for this to work in the e-free setting are given in (Al) and (A2) below. 

To continue, let (xo,/Uo), be fixed, their values to be determined later in (Al), and introdnce yo 
and zq by 


y = Co(x,z) +yo, 

z = y.o + zq. 

This leads to the following extended system obtained from (2): 

X = Fo{x,yo,zo) = /o(x) + Foy{x,yo, zo)yo + Fo^{x,yo, zo)zo, (11) 

yo = eo(x) -b ($(x) -b bo{x))zo + To(x, zo) -b (A(x) -b ao(x, zo))yo + Ro{x, yo, ^o), 
io = 0, 
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where 


Foy{x,yo,ZQ) = (using that F(x,Co(a;,/^o),^o) = 0) 

1 

dyF{x, Co(a;, ^o) + syo, + szQ)ds, 

1 

Foz{x,yo,ZQ) = I dzF{x,Co{x,y^o) + syo,fio + szo)ds, 

Jo 

eo(x) = G{x, Co(x, Mo), Mo), (12) 

$(x) = -dxCo{x, y,o)dzF{x, Co{x, Mo), Mo) + dzG{x, Co(x, Mo), Mo), (13) 

A(x) = -dxCodyF(x, Co(x, Mo), Mo) + dyG(x,Co(x,Mo),Mo), (14) 

and To = 0 (zq), oq = 0 ( 2 : 0 ), and Bq = 0(2/o)- Here we have just Taylor-expanded the right hands 
sides about (yo, zq) = (0,0). For later convenience we have also introduced /o = 0 and bo = 0. The 
subscripts on yo, zq and the functions fo, ■ ■ ■ ,bo, ■ ■ ■ , Ro are used to indicate that they later will be 
part of an iteration. The functions A and d* will not be updated. 

The function A = A(x) in (14) is precisely dyV ■ {—dxCo, 1) which according to (10) vanishes 
at a fold point (xo, Co( 2 : 0 , Mo)) at a given parameter value z = mo- Therefore we will assume that 
there is such a point where A vanishes. This only gives one condition on the pair (xo,Mo)- For 
the possibility of having a canard solution we need further conditions. These are contained in the 
following assumptions: 

(Al) The pair (xo,Mo) is so that A (I 4 ) and eo (12) vanish: 

-dxCo{xo, Mo)dyF{x, Co{xo, Mo), Mo)+dyG{xo, Co(x, Mo), Mo) = 0, 
G(xo,Co(xo,Mo),Mo) = 0. 



Given that A(xo) = 0 and eo(xo) = 0 we can define A and eo by 

A(x) = / dxA{xo + s{x — xo))ds, 

Jo 

respectively 

eo(x) = / dxBoixo + s{x - xo))ds, 

Jo 

so that A(x) = (x — xo)A(x) and eo(x) = (x — xo)eo(x). We further assume: 
(A2) The following non-degeneraey and “smallness” conditions hold true: Let 

<5o = \^o\uo, 

K^\A-X,. 

Then there exist e <C 1 so that 

60 < K-h, 

and 

I Foy \vQ, I Fqz \i/Q, I Tq \i,Q, I Rq I i,Q K e 


(15) 


(16) 


(17) 


Furthermore, either of the following conditions hold true: 
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(a) 

and 

(b) 

and 


mxo)\-^ K-h-\ 

Mxo)\-^ < K-h-\ 
I^Uo < K-h. 


(18) 

(19) 

( 20 ) 

( 21 ) 


Remark 2 We include case (b) in (A2) to cover the case where G = eH is small. 


The reason for supposing analyticity is that we in this case can present what we believe are 
optimal exponential estimates. However, we will not need any smoothness assumptions on how e 
enters beyond condition (Al) and (A2). In particular, the analysis of our method is not based on 
asymptotic expansions in e. 

The function cq in (15) is the obstacle to invariance of yo = 0, zq = 0; If cq = 0 then yo = 
0, zo = 0 corresponds to a canard solution. We have therefore introduced Sq in (16) as the error 
Sq = {eoliyQ = 0{e). The system (11) is our normal form. Our algorithm will be based upon 
applying transformations, affine in yo and zq, to (11) that seek to diminish the error cq. These 
transformations directly lead to a simple algorithm similar to that presented in [5] that we present 
in Algorithm 1. 

Remark 3 Suppose that G = eH is truly small and the slow and fast variables have been properly 
identified. Then there are known sufficient conditions for a canard explosion [3, 10]; 

(Bl) There exists a pair (xo,yo) so that 

dxCo{xo,flo) = 0, H{xoXo{xo,fio),lio) = 0. 


(B2) The following non-degeneracy conditions hold true; 

dyF{xo, Co(xo, Mo), Mo) / 0, d^F(xo, Co(xo, Mo), Mo) + 0, 
dzH{xo,C{xo,Mo),Mo) / 0,dxH{xo,Co{xo,Mo),Mo) / 0, 

Our condition (Al) replaces (Bl). They are equivalent when d‘^F{xo,Co{xo, Mo), Mo) / 0 by the 
implicit function theorem. Similarly, using (Al), it follows that the first three inequalities (B2) are 
equivalent to those in (A2), case (b), for e sufficiently small. In condition (A2), however, we do 
not require dxH{xo,Co{xo,Mo),Mo) / 0. This condition is included in (B2) because it guarantees 
that the nullclines of x and y are transverse at the fold point and that the equilibrium undergoes 
a Hopf bifurcation. It is the associated limit cycles that undergo rapid amplitude growth in the 
canard explosion. In agreement with [3], our algorithm and main result (Theorem 1 below) still 
apply without the need of this assumption but the results may have little dynamical significance. 
There is a well-known connection between the first Liapounov coefficient, the Hopf point and the 
canard point to lowest order in e [11, 3]. This was exploited in [12] as a numerical tool to estimate 
canard explosion points. □ 
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Remark 4 Suppose again that the slow and fast variables have been properly identified. Then one 
can actually replace A in (14) by —dxCodyF and ignore the small term dy{eH). Our main result 
still applies as the term ignored can be collected into oq. Indeed, the iterative lemma, Lemma 3, 
that is the basis of our proof of the main theorem, just assumes that ao(x, 0) is bounded from above 
by ce, see (41) below, for some c sufficiently large. However, in the general case, this term cannot 
be ignored. See also section 3.4 where we apply our algorithm to the Templator model. □ 

Remark 5 As described in Remark 4, the result of the paper still applies if |ao(x,0)|,yQ < |A| uo- 
Similarly, we can also allow for \fo\^Q <C |A|,^(, and |6 o|i,q <C |4>|,y(,. The iterative lemma. Lemma 3, 
still applies (see also (41) below). 

The modified iterative method for the computation of canard explosion. We are 

looking for a canard solution through a graph y = ^(x) for a value y of the parameter z. The 
invariance of the graph gives the equation (8) for which following Neishtadt [13] we wish 

to approach iteratively starting from y = Co(x,;Uo)) 2 : = yo and continue with solving the linear 
equations: 


Pn-l{x,z) + A{x)Cn{x,z) = 0, (22) 

Pn-i{x, z) = -dxC~^F{x, C~^,z) + G{x, C~^,z), (23) 

A(x) = -dxCodyF{xXo,l^o) + dyG{x, Co, Ao)- (24) 

Due to our assumption (Al), however, we are potentially dividing by zero as A(xo) = 0 when 


solving for the updates Cn- The idea of [5] is then to proceed by removing the singularity in the 
expression for (n by solving the equation pn-i{xo,z) = 0 (or z = y^ so that — pn-i{x, y^)/A{x) is 
well-defined. We collect this into an algorithm in the following: 

Algorithm 1 Suppose (Al) and (A2). To compute the canard explosion point do the following: 

1° Define A = A(x) to be 


A{x) = / dxA{xo +s{x-xo),yo)ds. 


Then by construction A(x) = (x — xo)A(x) where A is defined in (2f). 
2° Define 


eo{x)= dxPo{xo +s{x - xo),yo))ds. 
Jo 


Then by construction po{x,yo) = (x — xo)eo(x) where po is defined in (23)n=i- 

2° Iterate the following over n starting from n = 1 and = Co cind yfi = yo until \en\ has 
reached a desired tolerance: 

(i) Define 


Cn(x) — 


en-i{x) 

A(x) 


C^{x)=C^-\x) + Cn{x). 
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(a) Solve the following equation for z = fin-' 

Pnixo,n"'~^ + 2 :) = 0 , 

where 


Pn{x,l/" ^ + z) = -dxCF{x,C,h'' ^ + z) + G{x,C,h'' ^ + z). 


(Hi) Set 

pf^ = ^ + Pn, 

and let 

en{x)= dxPnixo +s{x - xo),p'^))ds. (25) 

Jo 

Then by construction pn{x,p^) = {x — xo)en{x). 

The graph y = C^{x) is then the approximation of the canard slow manifold, connecting repelling 
and attracting branches, at the explosion point z = p^. The error is described in Theorem 1. 

Main result. Our main result is contained in the following theorem which we prove in section 

4. 

Theorem 1 Suppose that the assumptions (Al) and (A2) hold true and that F and G are analytic 
in their arguments {x, y, z). (i) Fix first n > 0. Then, provided e is sufficiently small, the procedure 
defined in Algorithm 1 generates a sequence ofCi’s and pi’s so that 

n n 

satisfies (8) up to an error (25) of Cn = Moreover, (ii) there exists an N{e) = G 

N so that the error in (25) with n = N is exponentially small cn = 0{e~^F yjpp c > 0 and 
independent of e. 

The first part (i) only requires smoothness. The last part (ii) requires analyticity. 

Remark 6 The in Theorem 1 is in agreement with the results in [ 6 , 10] where the canard point 
is obtained as a smooth function of \/e. □ 


3 Applications 

3.1 Van der Pol 


In this section we consider the classical van der Pol system 

1 „ 

X = F{x,y) = y- -X + x, 
y = G{x, z) = e{z — x). 


This system has a canard explosion near z = pQ = 1 for e small. Asymptotic expansions yield a 
more accurate value 


p = l- 



173 3 
1024^ 




(26) 
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see e.g. [1, 17]. To use Algorithm 1 to compute the canard explosion point we first verify the 
conditions (Al) and (A2). Solving F{x,y) for y gives 

1 „ 

y = Co{x) = -X -X. 


Then 


A(x) = -da;CodyF(x, Co(x)) = 1-x^. 

We have a fold point at x = xq = 1 where 9xCo(a^o) = 0 where also A vanishes. To complete the 
verification of (Al) we must solve G(xo,z) = 0 for z = fiQ. We obtain /uq = 1. For (A 2 ) note that 

A(x) = -(l + x), (27) 


dividing x — xq = x — 1 out, so that A(xo) = —2 / 0 . Since <h = e we have to take case (b) in (A 2 ). 
The remaining assumptions can easily be verihed for e sufficiently small. 

We are now ready to apply Algorithm 1. For 1° we have (27) and for 2° we first note that 


Po(x,z) = e(z- x), 


and so eo = — e when dividing x — xq = x — 1 out from eo(x) = po(x, fio)- 

For 3° we first set Ci = —eo/^ = “^/(l + x), which finishes step (i) and define pi (ii) as 

Pi(x, z) = e(l + z) — e(l + x)“^(x^ + 3x^ + 3x + 1 — e). 

Setting pi(xo, z) = 0 gives z = pi = so that 

= 1 + = 1 - ^e, 


(28) 


correct to hrst order in e cf. (26). This finishes step (hi). We iterate this procedure and obtain the 
following approximations 


1 


27 




1 


173 


32 2048 

correct to order 2 respectively 3 in e. 


e3 + 0(e"), 


(29) 


3.2 Templator 

In this section we consider the Templator model [2, 4] 

X = F{x, y) = kuV^ + kry'^x - , (30) 

y = G{x, y,z) = z - kuy'^ - kry^x, 

and use Algorithm 1 to compute a canard explosion. Numerical computations indicate a canard 
explosion at 


p, = 0.419943, (31) 

[2, 4]. This system has no explicit e, yet numerical simulations show that the system exhibits a 
slow-fast structure. In [4] it is shown that various combinations of the parameters in the system 
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can locally be considered as small parameters, but no global parametrization in the form (2) exists. 
Here we proceed to find a canard explosion without any identification of an explicit small parameter. 
As in [4] we set ku = 0.01, /ct = 1, 9 = 1, and K = 0.05. We first solve F{x,y) = 0 for y = Co(2;) 
and obtain 

50v^2x(l + 150x + 5000x2) 

1 + 150x + 5000x2 




Note that (32) is independent of z. The equation also has a negative solution which we have 
discarded. We then realize that there is a point x = \/2/100 = 0.014142 where dxCo = 0- The 
function A(x) = —dxCodyF{x, (o) + dyG{x, Co, z) vanishes near this point aX x = xq = 0.014345. For 
this value of x = xq, we continue to verify the assumptions in (Al), and compute z = /iq giving 
G(xo, Co(xo), z) = 0. We obtain z = fio = 0.417681. The error is 0.5% in comparison with the value 
in (31). We then define A and cq by division of A and cq by x — xq- To verify (A2) we note that 

ch^l, 

A(xo) = 996.78, 
eo(a:o) = -17.157, 


so that 

= 0.017213. (33) 

A(xo) 

Since (33) is “small”, we are confident that (A2), case (a), is satisfied and we therefore proceed by 
applying Algorithm 1. Introducing Ci(3^) = —eo(x)/A(x) then gives the new error 

pi{x,z) = -dx(^F{x,C^) + G{x,C^,Po + z), C^ = Co + Cl- 

We solve pi(xo, z) = 0 for z = pi and obtain the improved approximation to the canard explosion 
point 

p,^ = po + Pi = 0.419883, 

an error of 0.01%. At the next step we get p'^ = 0.419938. The error is now 0.001%. 

3.3 Rotated van der Pol 

In this section we again consider the van der Pol equations, but we rotate the coordinates, replacing 
X by X — y and y by x + y: 

1 1 

X = F{x, y, z) = -(2x - (x - y f/3 + -e{p - (x - y)), 

y = G{x, y, z) = -^(2x - (x - yf/3) + ^e{p - (x - y)). 

mimicking a situation where the slow and fast variables have not been properly identified. We will 
demonstrate that the algorithm applies to this case too. The fold point in the original variables is 
no longer a fold point in the sense used in section 3.1 in the coordinates used here. Indeed solving 
for the x-nullcline gives 


y = Co(2;) = X — 6^/^x^/^ + 0(e) 


(34) 
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with non-zero derivative for the relevant x-values. The fold point now appears where A(x) = 
1 — vanishes. The 0(e)-term in A is following the discussion in Remark 4 not 

important: It can collected into oq = 0(e) in (11). We therefore re-define A as A(x) = 1 — 

Then A vanishes at the point x = 1/6, where y = —5/6-|-0(e) according to (34). The point 
(x,y) = (1/6,—5/6) is also the value obtained by transforming the fold point point (1,—2/3) in 
the original variables, used in section 3.1, to the current rotated variables. This alteration of A 
is in principle not needed: There is a point e-close to the point above where the old A vanishes. 
However, when we ignore this part then the calculations can actually be done by hand without 
using a computer algebra software. 

Note that 9a;A(xo) = —4 ^ 0 and we can therefore define 

32 

A(x) = -4 -b 4(x - xo) — :r{x - xq)^ + 0({x - xq)^). 

O 

Now we are in a position to apply Algorithm 1 2°. For this we first note that 

po(x, z) = -d^CoF(x, Co) + G(x, Co, z) = e6~‘^^^x~‘^/^(z - 6^/^x^/^), 
and we therefore obtain z = = 1 by solving po{xo, z) = 0. Following (hi), we then set 

eo(x) = -^e(3 - 18(x - xq) + 104(x - xq)^ + 0{(x - xq)^). 

O 

It is easy to verify the conditions in (A2). We therefore proceed as in Algorithm 1 (3°) starting from 

(1) setting Ci(3:) = —eo(x)/A(x) = —e6“^/^/(6^/^x^/^ + l)x^/^ which by construction is smooth also 
at X = Xq. We continue and obtain 

1 , 1 35 2 545 o 4, 

p -p, + p,-l--e--e e +0(e), 

2 1 3 2 18183 o 4, 

,..+^, = 1--,--, +0(t ), 

“ Ato + f,! + M2 + A,3 “ 1 - gf + 22*"^ “ ^ C(£‘*) 

correct to order 1, 2 respectively 3 cf. (26). Note that the corrections in (the terms — — 

|||e^) and p"^ (the term —are different from those in (28) and (29) being 0 and — 2518 ^^’ 
respectively, but the expressions are correct to the order expected by Theorem 1, even though the 
slow and fast variables have not been properly identified. 

3.4 Templator again 

As a final example we consider the Templator model (30) again with the parameters ku = 0.01, 
kr = I, q = I, and K = 0.05 as before. The solution of F(x, y) = 0 gave 

2 , , , 50W2x(l-bl50x-b 5000x2) 

« = =- l + 150x + 5000x^ -• 

having one single fold point at x = xq = \/2/100. There is, however, another canard explosion at 

p = 0.967555, (35) 

[2] that arises from another zero of A. Here we cannot ignore the term dyG in A as we did in 
the rotated van der Pol above. This term is not small as is illustrated in Fig. 1. Here it is made 
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visible that A has another zero at xq = 0.599393 where dx^ < 0. Therefore we can define A by 
A(x) = {x — xo)A(x). The zero of A gives rise to a singularity in the Eq. (22). To continue the 
verification of (Al) we note that solving po{xo, z) = 0 gives 

po = 0.967710. 

This gives a relative error of 0.01% in comparison with the value in (35). To verify the conditions 
(A2) we note the following: 

A(xo) = -3.6535, 
eo(xo) = -0.0521311, 

and <1> = 1 as above, so that 


= 0.014268. 

A(xo) 

This again gives us confidence to apply Algorithm 1. The expressions are quite messy so we leave 
out the details and just present the result of one iteration of 2° of the method: 

= 0.967560. 

The error is now 0.006%. A hnal additional application gives 

= 0.967558, 


reducing the relative error to 3 x 10 ® 



Figure 1: The graph of A (14) and A — dyG = —dxCodyF (dotted). 
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4 Proof of Theorem 1 


To prove our theorem we start from (11): 

X = Fo(x, 2/0, ^o) = fo(x) + Foy(x, yo, zo)yo + Fo^(x, yo, zo)zo, (36) 

yo = eo(x) + ($(x) + bo(x))zo + To(x, zq) + (A(x) + ao(x, zo))yo + Jio(x, yo, ^o), 
io = 0, 


where /o = 0 and bo = 0. From conditions (Al) and (A2), it is without loss of generality to take 
K = |A“^|jyp = 1. Indeed, we can introduce a new time r = K~^t to achieve this. Then also by 
(A2): 

4 = \e\uo < e < 1, (37) 

and furthermore, assuming case (a): 

i^(xo)r\ m < (38) 

Here (7$ > 0 is independent of e. The proof for case (b) where 

|^(xo)r^ <C$/e, 1^1 <C$e, (39) 

is almost identical (see Remark 7 below). 

We will seek to apply a sequence of transformations cpi, (/> 2 , ■ ■ ■■, (t>n to (11) that successively seek 
to diminish the Cj’s that appear as a result of these transformations. We define the transformations 
through an iterative lemma. 

The iterative lemma. To set up an iterative lemma we start from the normal form (36) with 
appropriate subscripts removed; 

X = F{x, y, z) = f{x) + Fy{x, y, z)y + F^{x, y, z)z, (40) 

y = e(x) + ($(x) + b{x))z + T{x, z) + (A(x) + a(x, z))y + i?(x, y, z), 

i = 0, 


where 


e = (x — xo)e(x), 

T = 0(z^) independent of y, and R = 0{y^). We suppose 

5 = \e\y, |a(x,0)|,,, |e“^6(xo)| < ce. 


(41) 


for some constant c > 0 and e sufficiently small. Note that these conditions are satisfied for our 
initial system (36). 

Then y = 0, z = 0 is an approximation to the canard solution. The accuracy of the approxima¬ 
tion is determined by 5. We will in the following apply a transformation (j) to (40), which will be 
based on two steps, seeking to improve the approximation of the canard solution. First, we define 
the solution, C, = C(^), of th® equation 

e(x) -|- A(x)C(x) = 0. (42) 


The solution is 


C(x) 


e(x) 


(43) 
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measuring 


Now, we set 


This gives 


\C\u<Kl 


y = C{x) + y+- 


X = f~^(x) + F^(x,y+,z)y+ + F+(x,y+,z)z, 
y+ = e+(x) + (4>(x) + b'^{x))z + r+(x, z) + (A(x) + a+(x, z))y + R'^{x, y, z), 
z = 0, 


with 


(44) 


f^{x) = f{x) + Fy{x, C, 0)C, 

F^{x,y+,z) = Fz{x,C + y+,z), 

e+(x) = -dxCf^ix) + a(x,0)C + i?+(x,C,0), 
b+{x) = b{x) - dxCFz{x,C,0) +dza{x,0)C + dzR+{x,C,0), 
a+(x,z) = a{x,z) - dxCFy{x,C,0) + dyR+{x,C,0), 

and where F+, T+ = 0{z‘^), which is independent of y, and R~^ = 0{y‘^) are determined by Taylor’s 
theorem. For e sufficiently small we then have 


(5+ = |e+ 




eC5 

~T 


for C sufficiently large, using (41), (44) and the fact that R'^{x,y,z) = 0{y‘^) is quadratic. 
The result is not yet appropriate for iteration as 


e+(xo) / 0. 


To account for this we transform 2 : by introducing z = y + z^ with y satisfying 

e^(a^o) + (^(a^o) + b^{xo))y + T^{xo,y) = 0. 


(45) 


By (38), and the contraction mapping theorem, there exists a solution 


y^ 


e+(xo) 

$(xo) ’ 


of (45), that satisfies 


\y\ < 2C^S+ < ( 46 ) 

Here we have used that = 0{z‘^), (41), the first estimate in (38) and the smallness of 5 and e. 
Then the resulting system reads 


X = F{x, C + y+,y + z+) = f+{x) + F+y{x, y+,z+)y+ + F+^(x, y+,z+)z+, 
y+ = e+{x) + {^{x) + b+{x))z+ +T+{x,z+) + (A(x) + a+{x, z+))y + R+{x,y, z+) 
z+ = 0, 


(47) 
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with 


f+{x) = f^ix) + F+ (x, 0, 

F+y{x,y+,z+) = F^{x,y+,fi + z+) 

e+(x) = e+(x) + (<h(x) + 6+(x))/U + T^{x,y,), 
b+{x) = 6+(x) +dzT+{x,fi), 


and where -F+z, T^, R+ are determined by Taylor’s theorem. By construction e+(xo) = 0. We 
estimate the new obstacle to invariance e+ of ?/+ = 0, = 0 as 


|e+|i/-^ < |e'’'|,y_^ + 2(7# 1^1 

^ eC(l + 2Cl)~d 


(48) 


using the second estimate in (38) and (41), for 5 and e sufficiently small. Since e+(xo) = 0 we write 
e+ as 


e+(x) = (x - xo)e+(x), e+(x) = / (9a;e+(xo + s(x - xo))(is, 

Jo 

and estimate our new error e+ on — 2^ using a Cauchy estimate and (48): 

eCil + 2Cl)6 






Remark 7 If we assume case (b) in (A2) then we obtain a €“^-factor in (46) so that |^| < 
2e“^|<I>(xo)|“^5^. However, we recover an estimate as in (48) using (39) (in place of (38) used 
above). □ 


We collect the results in the following iterative lemma: 

Lemma 3 Let ^ > 0 and = z/ — 2^ > 0. Then there exists a C > 0 so that the transformation 


(p : {x,y,z) (x+,y+,z+), 

X_|_ — 

y+ = C(x) - y, 

z+ = y- z, 

with C and y solving (43) and (45), respectively, maps (40) into (47) where 
|fl \b bjf-\i/_^_, Rjf- 11 /_|_, It T_|_I< C6, 

and 

6 + = le+|i.+ < ^ 6 , (49) 

provided e are suffieiently small. □ 
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The 0(e”) estimates in (i) of Theorem 1 follow directly from (49) as each application of the 
procedure introduces a factor of e. 

Exponential estimates. To obtain the exponential estimates we first apply the iterative 
lemma, Lemma 3, to (36) and obtain 






eCo 


So = 


16Co 


(50) 


using (37) and (49) with i/i = ^ setting here ^ = Co = ^ = ^1(1).^ Here Cq is the constant obtained 
from applying Lemma 3 to (36). Then we apply the lemma successively setting the measure of the 
domain reduction C in Lemma 3 to be^ 


c = V2eCoo, 


with 


Coo — 2 ( 70 . 

The (n + l)th application of Lemma 3 gives rise to a constant Cn- Then by (49) 

5n+l < < 2-^5n-l < 2-"5l, 

while Cn+i < Coo = 2(7o and 


(51) 


(52) 


l^n+l =l'n-2^ = -Uq - 2 nC > 0. 


(53) 


By the geometric sum formula and the fact that = O(e^), cf. (50), the requirement Cn < 2^70 
does not pose any restrictions on n for e sufficiently small. The only restriction on n is contained 
in the last inequality in (53); 


n < N{e) = 


1^0 



L4cJ 


_W2eCoo_ 


Here [uj denotes the integer part of u > 0. Hence Lemma 3 can be applied 0(e^/^)-many times, 
which by (52) with n = iV(e) gives the exponential estimate 


^N{e)+l < 2 


4v^2€Coo J ^ 


in Theorem 1. This then completes the proof. 
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